{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Interaction with other libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Keras\n",
    "\n",
    "- It's a very romantic notion to think that we can come up with the best features\n",
    "  to model our world.  That notion has now been dispelled.\n",
    "- Most *object detection/labeling/segmentation/classification* tasks now have\n",
    "  neural network equivalent algorithms that perform on-par with or better than\n",
    "  hand-crafted methods.\n",
    "- One library that gives Python users particularly easy access to deep learning is Keras: https://github.com/fchollet/keras/tree/master/examples (it works with both Theano and TensorFlow).\n",
    "- **At SciPy2017:** \"Fully Convolutional Networks for Image Segmentation\", Daniil Pakhomov, SciPy2017 (Friday 2:30pm)\n",
    "  - Particularly interesting, because such networks can be applied to images of any size\n",
    "  - ... and because Daniil is a scikit-image contributor ;)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Configurations\n",
    "\n",
    "From http://www.asimovinstitute.org/neural-network-zoo/:\n",
    "\n",
    "<img src=\"neuralnetworks.png\" style=\"width: 80%\"/>\n",
    "\n",
    "E.g., see how to fine tune a model on top of InceptionV3:\n",
    "\n",
    "<img src=\"inception_v3_architecture.png\"/>\n",
    "\n",
    "- https://keras.io/applications/#fine-tune-inceptionv3-on-a-new-set-of-classes\n",
    "\n",
    "\n",
    "- https://github.com/fchollet/keras/tree/master/examples\n",
    "- https://keras.io/scikit-learn-api/\n",
    "\n",
    "\n",
    "- In the Keras docs, you may read about `image_data_format`.  By default, this is `channels-last`, which is\n",
    "compatible with scikit-image's storage of `(row, cols, ch)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from keras.models import Sequential\n",
    "from keras.layers import Dense, Dropout\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "## Generate dummy data\n",
    "#X_train = np.random.random((1000, 2))\n",
    "#y_train = np.random.randint(2, size=(1000, 1))\n",
    "#X_test = np.random.random((100, 2))\n",
    "#y_test = np.random.randint(2, size=(100, 1))\n",
    "\n",
    "## Generate dummy data with some structure\n",
    "\n",
    "from sklearn import datasets\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "X, y = datasets.make_classification(n_features=2, n_samples=2000, n_redundant=0, n_informative=1,\n",
    "                                    n_clusters_per_class=1, random_state=42)\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.1, random_state=42)\n",
    "\n",
    "model = Sequential()\n",
    "model.add(Dense(64, input_dim=2, activation='relu'))\n",
    "model.add(Dropout(0.5))\n",
    "model.add(Dense(64, activation='relu'))\n",
    "model.add(Dropout(0.5))\n",
    "model.add(Dense(1, activation='sigmoid'))\n",
    "\n",
    "model.compile(loss='binary_crossentropy',\n",
    "              optimizer='rmsprop',\n",
    "              metrics=['accuracy'])\n",
    "\n",
    "model.fit(X_train, y_train,\n",
    "          epochs=20,\n",
    "          batch_size=128)\n",
    "score = model.evaluate(X_test, y_test, batch_size=128)\n",
    "\n",
    "print('\\n\\nAccuracy:', score[1]);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestClassifier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rf = RandomForestClassifier()\n",
    "rf.fit(X_train, y_train)\n",
    "rf.score(X_test, y_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f, (ax0, ax1, ax2) = plt.subplots(1, 3, figsize=(15, 5))\n",
    "\n",
    "mask = (y_train == 0)\n",
    "ax0.plot(X_train[mask, 0], X_train[mask, 1], 'b.')\n",
    "ax0.plot(X_train[~mask, 0], X_train[~mask, 1], 'r.')\n",
    "ax0.set_title('True Labels')\n",
    "\n",
    "y_nn = model.predict_classes(X_test).flatten()\n",
    "mask = (y_nn == 0)\n",
    "ax1.plot(X_test[mask, 0], X_test[mask, 1], 'b.')\n",
    "ax1.plot(X_test[~mask, 0], X_test[~mask, 1], 'r.')\n",
    "ax1.set_title('Labels by neural net')\n",
    "\n",
    "y_rf = rf.predict(X_test)\n",
    "mask = (y_rf == 0)\n",
    "ax2.plot(X_test[mask, 0], X_test[mask, 1], 'b.')\n",
    "ax2.plot(X_test[~mask, 0], X_test[~mask, 1], 'r.');\n",
    "ax2.set_title('Labels by random forest')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from keras.applications.inception_v3 import InceptionV3, preprocess_input, decode_predictions\n",
    "net = InceptionV3()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage import transform\n",
    "\n",
    "def inception_predict(image):\n",
    "    # Rescale image to 299x299, as required by InceptionV3\n",
    "    image_prep = transform.resize(image, (299, 299, 3), mode='reflect')\n",
    "    \n",
    "    # Scale image values to [-1, 1], as required by InceptionV3\n",
    "    image_prep = (img_as_float(image_prep) - 0.5) * 2\n",
    "    \n",
    "    predictions = decode_predictions(\n",
    "        net.predict(image_prep[None, ...])\n",
    "    )\n",
    "    \n",
    "    plt.imshow(image, cmap='gray')\n",
    "    \n",
    "    for pred in predictions[0]:\n",
    "        (n, klass, prob) = pred\n",
    "        print(f'{klass:>15} ({prob:.3f})')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from skimage import data, img_as_float\n",
    "inception_predict(data.chelsea())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "inception_predict(data.camera())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "inception_predict(data.coffee())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can fine-tune Inception to classify your own classes, as described at\n",
    "\n",
    "https://keras.io/applications/#fine-tune-inceptionv3-on-a-new-set-of-classes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SciPy: LowLevelCallable\n",
    "\n",
    "https://ilovesymposia.com/2017/03/12/scipys-new-lowlevelcallable-is-a-game-changer/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "image = np.random.random((512, 512))\n",
    "\n",
    "footprint = np.array([[0, 1, 0],\n",
    "                      [1, 1, 1],\n",
    "                      [0, 1, 0]], dtype=bool)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy import ndimage as ndi\n",
    "%timeit ndi.grey_erosion(image, footprint=footprint)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit ndi.generic_filter(image, np.min, footprint=footprint)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f'Slowdown is {825 / 2.85} times'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%load_ext Cython"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%%cython --name=test9\n",
    "\n",
    "from libc.stdint cimport intptr_t\n",
    "from numpy.math cimport INFINITY\n",
    "\n",
    "cdef api int erosion_kernel(double* input_arr_1d, intptr_t filter_size,\n",
    "                            double* return_value, void* user_data):\n",
    "    \n",
    "    cdef:\n",
    "        double[:] input_arr\n",
    "        ssize_t i\n",
    "        \n",
    "    return_value[0] = INFINITY\n",
    "    \n",
    "    for i in range(filter_size):\n",
    "        if input_arr_1d[i] < return_value[0]:\n",
    "            return_value[0] = input_arr_1d[i]\n",
    "    \n",
    "    return 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from scipy import LowLevelCallable, ndimage\n",
    "import sys\n",
    "\n",
    "def erosion_fast(image, footprint):\n",
    "    out = ndimage.generic_filter(\n",
    "            image,\n",
    "            LowLevelCallable.from_cython(sys.modules['test9'], name='erosion_kernel'),\n",
    "            footprint=footprint\n",
    "    )\n",
    "    return out"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.sum(\n",
    "    np.abs(\n",
    "        erosion_fast(image, footprint=footprint)\n",
    "        - ndi.generic_filter(image, np.min, footprint=footprint)\n",
    "    )\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit erosion_fast(image, footprint=footprint)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pip install numba"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Taken from Juan Nunez-Iglesias's blog post:\n",
    "# https://ilovesymposia.com/2017/03/12/scipys-new-lowlevelcallable-is-a-game-changer/\n",
    "\n",
    "import numba\n",
    "from numba import cfunc, carray\n",
    "from numba.types import intc, CPointer, float64, intp, voidptr\n",
    "from scipy import LowLevelCallable\n",
    "\n",
    "def jit_filter_function(filter_function):\n",
    "    jitted_function = numba.jit(filter_function, nopython=True)\n",
    "    \n",
    "    @cfunc(intc(CPointer(float64), intp, CPointer(float64), voidptr))\n",
    "    def wrapped(values_ptr, len_values, result, data):\n",
    "        values = carray(values_ptr, (len_values,), dtype=float64)\n",
    "        result[0] = jitted_function(values)\n",
    "        return 1\n",
    "    \n",
    "    return LowLevelCallable(wrapped.ctypes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "@jit_filter_function\n",
    "def fmin(values):\n",
    "    result = np.inf\n",
    "    for v in values:\n",
    "        if v < result:\n",
    "            result = v\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%timeit ndi.generic_filter(image, fmin, footprint=footprint)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parallel and batch processing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[Joblib](https://pythonhosted.org/joblib/) (developed by scikit-learn) is used for:\n",
    "\n",
    "\n",
    "1. transparent disk-caching of the output values and lazy re-evaluation (memoize pattern)\n",
    "2. easy simple parallel computing\n",
    "3. logging and tracing of the execution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.externals import joblib\n",
    "\n",
    "from joblib import Memory\n",
    "mem = Memory(cachedir='/tmp/joblib')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage import segmentation\n",
    "\n",
    "@mem.cache\n",
    "def cached_slic(image):\n",
    "    return segmentation.slic(image)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage import io\n",
    "large_image = io.imread('../images/Bells-Beach.jpg')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%time segmentation.slic(large_image)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%time cached_slic(large_image)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%time cached_slic(large_image)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[Dask](https://dask.pydata.org) is a parallel computing library.  It has two components:\n",
    "\n",
    "- Dynamic task scheduling optimized for computation. This is similar to Airflow, Luigi, Celery, or Make, but optimized for interactive computational workloads.\n",
    "- “Big Data” collections like parallel arrays, dataframes, and lists that extend common interfaces like NumPy, Pandas, or Python iterators to larger-than-memory or distributed environments. These parallel collections run on top of the dynamic task schedulers.\n",
    "- See Matt Rocklin's [blogpost](http://matthewrocklin.com/blog/work/2017/01/17/dask-images) for a more detailed example"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
